{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Concept\n",
    "==============\n",
    "\n",
    "\n",
    "$$\n",
    "\\newcommand{\\sumN}{\\sum_{n = 1}^N} \n",
    "\\newcommand{\\sumn}{\\sum_n} \n",
    "\\newcommand{\\bx}{\\mathbf{x}} \n",
    "\\newcommand{\\bbeta}{\\boldsymbol{\\beta}} \n",
    "\\newcommand{\\btheta}{\\boldsymbol{\\theta}} \n",
    "\\newcommand{\\bbetahat}{\\boldsymbol{\\hat{\\beta}}} \n",
    "\\newcommand{\\bthetahat}{\\boldsymbol{\\hat{\\theta}}} \n",
    "\\newcommand{\\dadb}[2]{\\frac{\\partial #1}{\\partial #2}} \n",
    "\\newcommand{\\by}{\\mathbf{y}} \n",
    "\\newcommand{\\bX}{\\mathbf{X}}\n",
    "\\newcommand{\\bphi}{\\boldsymbol{\\phi}}\n",
    "\\newcommand{\\bPhi}{\\boldsymbol{\\Phi}}\n",
    "$$\n",
    "\n",
    "## Model Structure\n",
    "\n",
    "*Linear regression* is a relatively simple method that is extremely widely-used. It is also a great stepping stone for more sophisticated methods, making it a natural algorithm to study first.\n",
    "\n",
    "In linear regression, the target variable $y$ is assumed to follow a linear function of one or more predictor variables, $x_1, \\dots, x_D$, plus some random error. Specifically, we assume the model for the $n^\\text{th}$ observation in our sample is of the form \n",
    "\n",
    "\n",
    "$$\n",
    "y_n = \\beta_0 + \\beta_1 x_{n1} + \\dots + \\beta_Dx_{nD} + \\epsilon_n. \n",
    "$$\n",
    "\n",
    "\n",
    "Here $\\beta_0$ is the intercept term, $\\beta_1$ through $\\beta_D$ are the coefficients on our feature variables, and $\\epsilon$ is an error term that represents the difference between the true $y$ value and the linear function of the predictors. Note that the terms with an $n$ in the subscript differ between observations while the terms without (namely the $\\beta\\text{s}$) do not."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The math behind linear regression often becomes easier when we use vectors to represent our predictors and coefficients. Let's define $\\bx_n$ and $\\bbeta$ as follows:\n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "\\bx_n &= \\begin{pmatrix} 1 & x_{n1} & \\dots & x_{nD} \\end{pmatrix}^\\top \\\\\n",
    "\\bbeta &= \\begin{pmatrix} \\beta_0 & \\beta_1 & \\dots & \\beta_D \\end{pmatrix}^\\top.\n",
    "\\end{align}\n",
    "$$\n",
    "\n",
    "Note that $\\bx_n$ includes a leading 1, corresponding to the intercept term $\\beta_0$. Using these definitions, we can equivalently express $y_n$ as \n",
    "\n",
    "\n",
    "$$\n",
    "y_n = \\bbeta^\\top \\bx_n + \\epsilon_n.\n",
    "$$\n",
    "\n",
    "\n",
    "Below is an example of a dataset designed for linear regression. The input variable is generated randomly and the target variable is generated as a linear combination of that input variable plus an error term."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 327,
   "metadata": {
    "tags": [
     "hide-input"
    ]
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYEAAAEeCAYAAABsaamyAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3deXxU1d3H8c8vBAKhuBQFLZQqLoALSam7UkRpBVu1faxFrCgVagURN8S1rUi1oKhoBRW3VgSpOy4VIY8oi1q2B1A2iVZZfB62AgIh+3n+OEkMmMAkmTt3Zu73/XrlhZm5M/dkTO7v3nPP+R5zziEiItGUEXYDREQkPCoCIiIRpiIgIhJhKgIiIhGmIiAiEmEqAiIiEaYiICISYSoCIiIRpiIgsg9mNtXMnJn9OcB9/MrMXjazL81sl5mtNLO/mFmLoPYpAioCIntlZn2AnATsaihQBtwG9AQeBQYC081Mf6cSGP1yidTCzA4AHgRuqOfrvzCzO2Pc/Dzn3K+dcxOdc+8758YAQ4CTgTPrs/84tk3SmIqAxJWZTTezD2t4/HgzKzGzS8JoVz3dCyx1zj0f9I6ccxtreHhexb9tqj9oZkdWfJbD93j8UTPbbmYnBNFGM2tvZn8zs3VmVmxmG8xshpk1DmJ/khiZYTdA0s5s4DYzy3LOFQGYmQHjgA+cc5Oqb1zxXKMY3tc558ri3tpamNkZwGUkpiuoNt0q/l1e/UHnXL6ZPQlcb2Z/dc5tMrM/AlcAP3POzY93QyquiuYA/8J3U20Bvgu0cc6VxHt/kjgqAhJvc4AmwA+Bjyoeuww4BehSw/bdgBkxvO/7BNAtUpOKM9vHgdHOuZUxvqa2YpZhZtX/zmIqZmbWBrgLyKvloD4c/7nebGYrgD8BfZxzeQG17TTgEOAfQB5Q7Jwr3dfPIcnPFCUt8WRm3wG2AkOdc2MqziBXApOcc9fXsH0LoEMMb709lgOymfUApsfwfu87586s5T3uwJ9VH+uc21XxmAPuds7dUctrziTGYlbbfqu913eA94DvASc559bWst3dwI34k7lrnXNjg2qbmR2ML/BHVTz0H+dcyxjeU5KcrgQkrpxzO8xsMf7MH+BuoBx/plqTHcCiWN46xiZ8AHSKYbuCmh40s3bA7cAAIMvMsqo9nVVR1LbXcMa8ADhxj8deB94Exld7bPveGmVmTSte1x7oVlsBqLAKyAJm11YA4ti25sDbwAPAEnyhl3TgnNOXvuL6BTwMfIHv/ikDLt3LtmfiD/D7+novQW2PpT25Mb7XF8Cdddh3Y+AtfGE8ZR/bngUU4YteOZBTx58z5rYB+1Vsf17Yv1v6iv+XrgQkCLOBa4BngTnOuef2sm1NZ6k12esZdBwtArrX8PgM4DngKSA/3jutmAswETgbf3P3o71s2wV4DXgSuB74FLgH+Fm821Xhx8AP8N16kmZUBCQIcyr+7UjNN4OrOOe2A3EfzVJfzrmt+P743fh7q3zpnPvWc3EyFrgI332208xOqfbcWlfRLWRmR+K7ZaYB1zjnyiuGij5tZj92zs0MoG2Vw1efNbP7gM1AW+As59wVAexPEkjzBCQIO4Bi4BHn3JKwG5MielX8ezvw4R5fAwDM7BD8wX858BvnXHnFa54FVgAjg2iYc+5fQF98d9UzwD/xM5s/D2J/klgaHSRxZ2b3A5cAHZ1z28Juj4jUTt1BEhdmlo2fWNUVuBa4SAVAJPmpCEi89ACmAOvwY9ZfDbk9IhIDdQeJiESYbgyLiERYSnUH9ezZ002dOjXsZoiIpBqr7YmUuhLYtGlT2E0QEUkrKVUEREQkvlQEREQiTEVARCTCVARERCJMRUBEJMJUBEREIiyl5gmIiERFQVEp24tK2VJQzIHZTWiRlUl2VvwP2SoCIiJJpqColLzl6xn64hKKy8pp0iiD0Rd1pken1nEvBOoOEhFJMtuLSqsKAEBxWTk3vbSE7UWlcd+XioCISJLZUlBcVQAqFZWWs7WgJO77UhEQEUkyB2Y3oUmj3Q/PWZkZHJDdOO77UhEQEUkyLbIyGX1RZ7Iy/SE6KzOD+36VQ4sAbgyn1HoCJ5xwgps/P2nWJBcRCUzl6KCtBSUckN24oaODak0R1eggEZEklF1x0G+9X9NA96PuIBGRCFMREBGJMBUBEZEIUxEQEYkw3RgWEUlSicgPUhEQEUlCicoPUneQiEgSqswParX5K45d/1lg+UEqAiIiSWjrlu1cOWsSeU8N4u53xoJzgeQHqQiIiCSbvDyO6HE6Q2c9R94RJ3HVL24Ds0Dyg3RPQEQkWXz1FdxwA/zjH2QccSQfjn2OG7/6LkWl5YHlB6kIiIiErbQU/vpX+NOfoLgYhg8nY9gwciyTmfHLD6qRioCISJjmzIFBg2DJEujVyxeDI44AIBsCzw/SPQERkTBs3Aj9+8MZZ8CWLfDKK/DWW1UFIFFCLQJm1tPMVppZvpndEmZbREQSorwcxo+HDh3g2Wfh5pth+XL45S/Bak18DkxoRcDMGgFjgV7AMUAfMzsmrPaIiARu4UI47TT4/e+hc2dYvBhGjoTmzUNrUphXAicB+c65z51zxcBk4IIQ2yMiEoxt22DIEDjxRPj3v2HCBJgxA44J/7w3zCLQBlhT7fu1FY/txsyuNLP5ZjZ/48aNCWuciEiDOQcTJ/qun7FjYeBAWLkSLr00lK6fmoRZBGr6BL611qVzbrxz7gTn3AkHH3xwApolIhIHy5fD2Wf7A367djB3LjzyCBxwQNgt202YRWAt8P1q37cFvgqpLSIi8bFzJ9x6K+TkwKJF8Nhj8OGH8KMfhd2yGoU5T2AecJSZHQ6sAy4GLgmxPSIi9eccvP667/tfvRr69YNRo6BVq7BbtlehFQHnXKmZDQbeARoBTzvnlobVHhGRevv3v/3B/8034bjjYNYsP/4/BYQ6Y9g590/gn2G2QUSk3oqK4L774O67ITMT7r8frrkGGsc35C1Iio0QEamP6dNh8GD49FO46CJ44AFo2zbsVtWZYiNEJFAFRaWs/7qQFf/3Neu/LqQgzouiJNy6dXDxxfDTn/rZv1OnwgsvpGQBAF0JiEiAErVEYkJUJn3+8Y9QUgLDh8OwYdA0uHC3RNCVgIgEpnKJxOKycoDAlkgM3Jw50KWLz/r/8Y9h6VJfDFK8AICKgIgEaEtBcVUBqBTEEomB2bgRrrjCj/TZutUnfb75ZsKTPoOkIiAigTkwuwlNGu1+mAliicS4Ky+Hxx/3cQ8TJvhunxCTPoOkIiAigWmRlcnoizqTlekPNUEtkRhXCxfCqafCVVf5Wb+LF/tJXyEmfQYpif9PiEiqy87KpEen1swc1j3QJRLjYutW+MMfYNw4OPhgfwXwm9+k3Zn/npLw/4SIpJPsioN+kEskNohzMGkS3HijvwcwaBCMGJF0QW9BUREQkehavtwf9N97D046yS/vmKRBb0HRPQERiZ7KpM/K1b2SPOkzSLoSEJHocA6mTIFrr/VJn7/9rb/pG+G1SlQERCQaPv/cJ32+9VbKJX0GSd1BIpLeior8jd5jj4X33/dJnwsXqgBU0JWAiKSv6dPh6qth1aqUTvoMkq4ERCT9rFsHvXv7pE/nUj7pM0gqAiKSPkpK/Nl+x45+qcfhw+Hjj+Gcc8JuWdJSd5CIpIfZs/2Y/48/hnPP9bHP7dvX6S0KikrZXlTKloJiDsxukryzm+MovX86EUl/GzfCzTfDM8/A978Pr74KF1xQ57iHtFr7oA7UHSQiqam2pM9f/KJeeT9ps/ZBHakIiEjqWbAg7kmfKb/2QT2pCIhI6ti6Fa65xuf8fPklPPccvPsuHHNMg986Zdc+aCAVARFJfs75A37Hjj7qedAgWLEirlHPKbn2QRyk908nIqlv2TLKrhpIo1kz2fXDH1Hwwqs0O/nEuN+sTam1D+JIVwIikpx27oRbbsHl5FC2eDF/6HUNx/zkT5z6zlbylq+nIIAbtpXrHnQ4pAWt92ua9gUAdCUgIsmmMulzyBBYs4bCSy/jzJa9WN+0BfDNqJ2Zw7pH4iAdNF0JiEjy+PxzOO88v6D7AQfA7Nl8ed9fqwpApSiM2kkUFQERCV9NSZ8LFsDpp0d21E6iqAiISLimTYPjj4c//hHOP9+P+rnhBmjsD/JRHbWTKPoURSQc69b5g/0LL8BRR8E77/jUzz1EddROooTyKZrZfcB5QDHwGfBb59zWMNoiIglWUuLD3f70Jygt9d1AN90EWVm1viS74qDfer+mCWxoNITVHTQdOM451xn4FLg1pHaISCLNnu0Xc7/xRujWDZYuhTvu2GsBkGCFUgScc9Occ5WDfD8CtNKDSDrbuNEv6t61K2zbBq+9Bm+8UeeoZ4m/ZLgxfAXwdtiNEJEAlJV9k/T53HNwyy2wbFm9op4lGIHdEzCzPOCQGp663Tk3pWKb24FSYOJe3udK4EqAdu3aBdBSEQnEggU+42fuXDjzTBg7Ni5BbxJf5pwLZ8dmlwNXAWc75wpiec0JJ5zg5s+fH2zDRKRhtm71/fyPPgoHH+yXe+zTR2f+4ar1ww9rdFBP4GagW6wFQESSnHMwcaK/6btpk78KGDHCz/yVpBXWQNtHgCxguvmzg4+cc1eF1BYRaahly+Dqq+G993zW/9tvQ5cuYbdKYhBKEXDOHRnGfkUkznbu9Gf7998PLVr4m8ADBkBGMow5kVhoyp2I1N0eSZ9ccQWMHOnvAUhKUbkWkbqpIemTp55SAUhRKgIiEps9kz4feAAWLoTTTw+7ZdIA6g4SkX2bNg0GD4ZVq6B3b38PoE2bsFslcaArARGp3bp18Otfwznn+O/feQcmT1YBSCMqAiLybSUlvrunY0ef8XPXXfDxxzVGPUtqU3eQiOxu9mw/0evjj+Hcc33ss4Le0pauBETEqynp8803U74AFBSVsv7rQlb839es/7qQgqLSfb8oQnQlIBJ1ZWXwxBNw222wfbtP+rzjDmjePOyWNVhBUSl5y9cz9MUlFJeV06RRBqMv6kyPTq21MlkFXQmIRNmCBXDqqTBwIOTkwOLF8Je/xLUAhHkmvr2otKoAABSXlXPTS0vYrquBKiqFIimuoKiU7UWlbCko5sDsJrGtv1uZ9DluHLRqBRMmwG9+E/ekz7DPxLcUFFcVgEpFpeVsLSjRUpUVdCUgksIqD7JdR82g55hZdB01g7zl62s/23bOL+7SoYOPeh48GFasgEsvDSTqOewz8QOzm9Ck0e6HuazMDA7IbpyQ/acCFQGRFFang+zSpdC9O/TtC4cdBvPmwcMPBxr1vLcz8URokZXJ6Is6k5XpD3VZmRnc96scWuh+QBV9EiIpLKbujh07fNzDAw/4pM/x46F//4QkfVaeiVdvYyLPxLOzMunRqTUzh3Vna0EJB2Q3jq27LEJ0JSCSwvba3eEcvPqqX9Lx3nvhsstg5Ur43e8SFvWcDGfi2VmZtN6vKR0OaUHr/ZqqAOwhtOUl60PLS4rsrvKewE0vLaGotLzqIPuTJttpNvR6v7jL8cf7/v+Qgt4qb1zrTDxUtd7wUREQSXG7HWQzyjjwkTE0uXckNG7s4x6uuQYyddCNuORaY1hE4ie74sy69Ucz/RKP+fk+9O2BBxT0JvukewIiqW7t2m+SPs187PM//qECIDFRERBJVZVJn506+aTPESN86NtPfhJ2yySFqDtIJBVVT/r82c980ufhh4fdKklBuhIQSSU1JX2+8YYKgNTbPouAmV1mZpvNLGuPxyea2evBNU1EqpSVweOP+7iHiRN90ueyZXDBBYHEPUh0xHIl8GLFdhdUPmBm+wO/BJ4KqF0iUmnBAjjtNLjqqsCSPiW69lkEnHO7gInAFdUevgT4GngroHaJyNatPuDtpJPgyy/9FcC77/obwSJxEus9gSeAn5hZ24rvrwD+7pxTKLdIvDnno50rkz4HDfJJn5dcoq4fibuYRgc55xab2UKgn5m9BpwAXBpoy0SiaOlSf9CfORNOPtnHPnTpEnarJI3VZXTQE0A/YAAwxzm3MpAWiUTRjh1w882QmwuffOKTPj/4QAVAAleXIvA8cAgwEN0QFokP5+CVV3w//733wuWXJzzpU6It5t8y59x24AWguOJfEWmIzz6Dn/8cLrwQvvtdPwHsySfhoIPCbplESF1PNQ4FJjvndgbRGJFIKCz06Z7HHuv7/h980A8DDSnqWaItphvDZvZdoAfwUyAnng0ws6HAfcDBzrlN8XxvkerqtSB7vL3zjh/2mZ8PvXvD/fcr6E1CFetfwELgu8BtzrlP4rVzM/s+8BNgdbzeU6QmlYuvVK7H26RRBqMv6kyPTq0TUwjWroXrr4eXXoKjjvJJnwp6kyQQU3eQc+4w59x+zrlRcd7/g8AwIHVWtpGUVKcF2eOppMSf7XfqBG++CX/+s5I+JamEliJqZucD6yrmIOxtuyuBKwHatWuXoNZJuolpQfZ4mz0bBg70Qz6V9ClJKtAiYGZ5+GGle7oduA1/j2GvnHPjgfHgl5eMawMlMioXZK9eCKoWZI+3DRv8mP+//Q3atfNJn+efr9m+kpQCHYjsnOvhnDtuzy/gc+BwYLGZfQG0BRaaWU0FQ6TBWmRlMvqizmRl+l/5ygXZW8TzfkBZGTz2mJI+JaUkxULzFYXghH2NDtJC89IQuy3Int04vqODFizwXT/z5kH37jB2rILeJJlooXmRqgXZ43kPYOtWuOMOGDcOWrWC555T0JuklKQoAs65w8Jug0idOOcP+EOHwqZNfuz/iBGw//6B7C4p5jhIWtJvkUhdJTjpM/Q5DpLWlFAlEquQkj5Dm+MgkaAiILIvISd97m2Og0hDqQiI7E1+vp/oVZn0OWdOwpM+K+c4VBfYHAeJHBUBkZoUFsLw4XDccTBr1jdJn6edlvCmJGSOg0RWUswTiJXmCUhCTJ3qR/t89lnSJH0GOsdBokDzBET2qXrS59FHw/Tp0KNH2K0CAprjIIK6g0R80ufo0dCxo0/6HDEClixJmgIgEiRdCUi0zZrl4x6WLvVLPT78sJI+JVJ0JSDRtGED9OsHP/6xH/8/ZQq88YYKgESOioBES/Wkz0mT4NZb/VXA+eeH3TKRUKg7SKJDSZ8i36IrAUl/W7fC1VfjTjyRsi9X89UjT7D+1bcoaH9U2C0TCZ2KgKQv52DCBOjQAffYY/y7dz9OvPQRTltzKF3vfY+85espUP6ORJyKgKSnpUvhzDPhssvg8MP5z3tz6Nn+Iv7TuBmgEDaRSioCkl527IBhw76V9LnxqGMUwiZSA90YlvRQmfR53XV+5m///jByZFXQW0IXmhdJIboSkNT32Wdw7rnwq19By5Y1Jn0qhE2kZgqQk9RVWAijRsFf/gJNmvi4h6uvhsyaD+wKYZMIU4CcpJnqSZ8XX+yTPr/3vb2+RCFsIt+m7iBpsIKiUtZ/XciK//ua9V8XBjvscs0a3+3Tqxc0auSTPp9/fp8FQERqpisBaZCELYJeUgIPPQR33umjH+6+G268EbKy4rcPkQjSlYA0SEIWQZ81C374Q7jpJh/3sGwZ3HabCoBIHKgISIMEugj6HkmfRS+9wvqJL7KiWcvgu51EIkLdQdIggYy/Lyvzk7xuuw127oRbb6Vg6M3krd7B0FEzgu12EokYXQlIg8R9/P2CBXDqqTBoEHTp4lf4uucetmdmBd/tJBJBOoWSBsnOyqRHp9bMHNa9YePvt2yBO+6ARx+F1q1h4kTo0wfMD2/eW7eThnyK1J+KgDRYg8bfOwfPPQdDh8KmTXDNNXDXXbD//rttptgHkWCoOygNJXTcfkPskfTJ/Pl+GOgeBQAU+yASFMVGpJmEjdtviB07/Nn+gw/Cfvv56IcrroCMvZ+TKPZBpN5qjY3QlUCaSci4/fpyDl5+2S/peN99cPnlsHIlDBiwzwIAVHU5dTikBa33a6oCIBIHoRUBM7vGzFaa2VIzuzesdqSbQMftN0R+/j6TPusiZbq8RJJcKKdSZtYduADo7JwrMrNWYbQjHSXdDdQ9kz4fesgP/6wl6TMWKdHlJZIiwroSGAiMdM4VATjnNoTUjrSTVDdQp06F447zeT+//CWsWAFDhjSoAECSd3mJpJiwTpuOBrqa2d1AITDUOTevpg3N7ErgSoB27dolroUpKm7j9htizRq4/nrf/9+hA+Tlwdlnx+3tNWdAJH4COzKYWR5wSA1P3V6x3wOBU4ATgRfMrL2rYaiSc248MB786KCg2ptOQsvNLymBMWNg+HAoL4c//9mP/49z0FvSdXmJpLDAioBzrkdtz5nZQOCVioP+XDMrBw4CNgbVHgnYzJm+r3/pUjjvPHj4YTjssEB2VdnlddNLSygqLdecAZEGCOuv5jXgLOA9MzsaaAJsCqkt0hAbNviI52efhR/8AKZMgfPPD3SXDenyqpxrsKWgmAOzm2iugUReWL/9TwNPm9knQDFweU1dQZLEysrg8cfh9tt90udtt/n/zs5OyO7r0+WlUUUi3xbKb75zrhi4NIx9SxzMmwcDB/rEz7POgrFjoWPHUJsUyxl+baOKZg7rriIgkaXffIndli3+bP+xx3zS5/PPQ+/eVUmfYYn1DF+jikS+TbERsm/O+T7/Dh18F9CQIX7M/8UXh14AIPZ5A5WjiqrTqCKJOhUB2btPPoFu3XzOzxFH+C6gMWNqTPoMS6xRGUk1kU4kSei3X2q2Z9LnE0/ElPRZF/EaqRPrvIGkmEgnkmT02y+7cw5eeQWuuw7WroX+/WHkyHoHvdUmniN16jJvILSJdCJJSusJyDfy8/3KXlOnQk4OjBsHp50WyK7Wf11I14pF4ytlZWYwc1j3eh2gtdaAyF7VevNOfyXy7aTPMWPg6qsbHPS2N/EeqaMzfJH6URGIuqlTYfBg+Owzv7D76NHwve8Fvlvl/4gkB40Oiqo1a/wCL716+TP+vDyYNCkhBQA0UkckWeieQNTsmfR5xx1w441xT/qMhfrxRRJG9wSEhCZ9xkL9+CLhU3dQFKxf7yd7devmx/9PmQKvvx5qAUgnWu9YUpmuBNJZWRmMH+8TPnfuhFtv9d0/CUr6jAIlk0qq05VAupo3D045xXf//OhHsGQJ3HOPCkCcab1jSXUqAulmyxZ/4D/5ZD/jd9IkmD499KjndBVrbpFIslIRSBe1JX326ZMUSZ/pSsmkkupUBNLB0qVJn/SZrjTfQVKd5gmksh07/Hj/MWN80ue998JvfxvXpE/ZN813kBSgeQJpZc+kzwEDfO5PnJM+JTaa7yCpTKeMqSY/H84910c+tGwJH3zgs/5VAESkHlQEUkVhIdx5Jxx3HMyZAw89BPPnw6mnht0yEUlh6g5KBW+/7XP+P/vMr+t7//0JC3oTkfSmK4FktmYNXHih7/6pTPp8/nkVABGJG10JJKM9kz7vvju0pE+RhiopKWHt2rUUFhaG3ZS017RpU9q2bUvjxrHPU1ERSDbVkz7PP9/3/SvoTVLY2rVradGiBYcddhimiYuBcc6xefNm1q5dy+GHHx7z69QdlCzWr4fLLvOTvnbu9CmfU6aoAEjKKywspGXLlioAATMzWrZsWecrLhWBsJWV+QXdO3SAyZPh9tu/yfsXSRMqAIlRn89Z3UFhmjcPBg70MQ9nnw2PPKKgNxFJKF0JhKF60udXX/kRP0r6FAnE5s2byc3NJTc3l0MOOYQ2bdpUfV9cXBy3/eTl5WFm/P3vf696bN68eZgZY8aMifl98vPzyc3NbfA2sdKVQCJVJn3edBNs3uyTPu+6y+f+iEggWrZsyaJFiwC48847+c53vsPQoUN328Y5h3OOjAbmbh1//PFMnjyZyy+/HIDJkyeTk5PToPcMWiSKQGXA15aCYg7MbhJOwNcnn/iz/1mz/CzfadMgTpVcJGVcdx1UHJDjJjfXD6muo/z8fH7xi19wxhln8K9//YvXXnuNnJwctm7dCvgDeF5eHk8++STr169n4MCBrF69moyMDB5++GFOOeWUb71n+/bt2bhxI5s2baJly5ZMnz6dXr16VT2/cOFCBg4cyK5duzjqqKN4+umn2X///Zk3bx79+/enefPmnH766VXbl5aWMmzYMGbPnk1hYSFDhgxhwIAB9fiQahdad5CZ5ZrZR2a2yMzmm9lJQeyncvm/rqNm0HPMLLqOmkHe8vWJWwd2xw5/5p+b62/4PvkkzJ6tAiCSBJYtW0b//v35n//5H9q0aVPrdkOGDGHYsGHMnz+fF154Ya8H4gsvvJCXXnqJmTNncvLJJ+82Zv/SSy/l/vvvZ8mSJXTo0IERI0YA0K9fPx599FE+/PBDysrKqrYfP348rVq1Yu7cucybN4+xY8eyevXqOPzk3wjzSuBeYLhz7m0zO7fi+zPjvZPalv+bOax7sFcDlUmf114L69Yp6VME6nXGHqQjjjiCE088cZ/b5eXlsXLlyqrvt2zZwq5du2jWrNm3tu3duzd9+/bl6KOPpk+fPrz77ruAvzdRWFjIGWecAcDll19O37592bRpE7t27aq6Aujbty8zZswAYNq0aSxfvpzJkycDsG3bNlatWsUPfvCDhv3g1YRZBBxQ2Rm+P/BVEDvZ2/J/gUX/5uf7rJ+pUyEnB158UUFvIkmoefPmVf+dkZFB9fVVqo+3d84xd+5cmjRpss/3bNOmDc453n//fcaNG1dVBPa2dkttQzudc4wbN46zzz57t8fz8/P32Y5YhTk66DrgPjNbA4wGbq1pIzO7sqK7aP7GjRvrvJOELv+3a5eSPkVSVEZGBgceeCCrVq2ivLycV199teq5Hj16MHbs2KrvF+3jvsaIESMYNWrUbjeaDzroIJo1a8YHH3wAwIQJE+jWrRsHHXQQTZs25cMPPwRg4sSJVa8555xzGDduHKWlvvt65cqV7Nq1q+E/bDWBXgmYWR5wSA1P3Q6cDVzvnHvZzH4NPAX02HND59x4YDz4lcXq2obK5f9uemkJRaXlwS3/Vz3ps08fn/R56KHx3YeIBGrUqFH07NmTdu3accwxx1BUVATA2LFjGThwIM888wylpaV07959t6Kwp8ounz1NmDCh6sbwkUceyTPPPAPAM888w4ABA83pbhcAAAdBSURBVGjevDk//elPq7b//e9/z+rVq6uGg7Zq1YopU6bE68cFQlxe0sy2AQc455z5a6Ftzrm9jpWs7/KSgS7/t2aNH/Hwyit+1u/YsX7il4gAsHz5cjp16hR2MyKjls87KZeX/AroBrwHnAWsCmpHgSz/V1z8TdKnc3DPPXDDDUr6FJGUEmYR+B3wkJllAoXAlSG2pW7ef9+P+V+2TEmfIpLSQisCzrnZwI/C2n+9rF/vx/xPmOAP+q+/rqA3EUlpyg6KhZI+RSRNRSI2okHmzvVdPwsWQI8ePumzQ4ewWyUiEhe6EqjNli0+5vmUU75J+pw2TQVARNKKisCenIO//c0f7MeP97EPK1bAxReDFsYQSUmNGjUiNzeXY489lpycHB544AHKy8v3+povvviCSZMmJaiF4VF3UHWffOLP/mfPVtKnSEiCSP1t1qxZ1SzfDRs2cMkll7Bt2zaGDx9e62sqi8All1zSoH0nO10JAGzfDkOH+gP+8uXw1FNK+hQJQSJSf1u1asX48eN55JFHcM7xxRdf0LVrV7p06UKXLl2qYh1uueUWZs2aRW5uLg8++GCt26W6aF8JOAcvv+xn/K5bB7/7nU/6bNky7JaJRFKiUn/bt29PeXk5GzZsoFWrVkyfPp2mTZuyatUq+vTpw/z58xk5ciSjR4/mzTffBKCgoKDG7VJddItAfj4MHgzvvOPP+F96yd8EFpHQJDL1tzIyp6SkhMGDB7No0SIaNWrEp59+WuP2sW6XaqJXBHbtglGjYORIH/Hw0EN+CGhm9D4KkWRTmfpbvRAEkfr7+eef06hRI1q1asXw4cNp3bo1ixcvpry8nKZNay42Dz74YEzbpZpo3RP45z99zPPw4fBf/+VH/QwZogIgkiQqU3+zMv2hKYjU340bN3LVVVcxePBgzIxt27Zx6KGHkpGRwYQJE6pW9mrRogXbt2+vel1t26W6aBz9nIPevf3iLh07wn//N5x1VtitEpE9ZGdl0qNTa2YO6x7X1N9du3aRm5tLSUkJmZmZ9O3blxtuuAGAQYMGceGFF/Liiy/SvXv3qoVmOnfuTGZmJjk5OfTr16/W7VJdaFHS9VHfKGkA/vAHaN7cJ33GsDqQiMSHoqQTK5WipBOrYkFnERH5RrTuCYiIyG5UBEQkcKnU7ZzK6vM5qwiISKCaNm3K5s2bVQgC5pxj8+bNdR66Gp17AiISirZt27J27Vo2btwYdlPSXtOmTWnbtm2dXqMiICKBaty4MYcffnjYzZBaqDtIRCTCVARERCJMRUBEJMJSasawmW0Evgy7HXVwELAp7EaETJ+Bp89Bn0GlMD6HTc65njU9kVJFINWY2Xzn3AlhtyNM+gw8fQ76DCol2+eg7iARkQhTERARiTAVgWCND7sBSUCfgafPQZ9BpaT6HHRPQEQkwnQlICISYSoCIiIRpiIQMDO7z8xWmNkSM3vVzA4Iu02JZmYXmdlSMys3s6QZGpcIZtbTzFaaWb6Z3RJ2e8JgZk+b2QYz+yTstoTFzL5vZjPMbHnF38K1YbepkopA8KYDxznnOgOfAreG3J4wfAL8FzAz7IYkkpk1AsYCvYBjgD5mdky4rQrF34AaJypFSClwo3OuE3AKcHWy/C6oCATMOTfNOVda8e1HQN1yXtOAc265c25l2O0IwUlAvnPuc+dcMTAZuCDkNiWcc24m8J+w2xEm59z/OucWVvz3dmA50CbcVnkqAol1BfB22I2QhGkDrKn2/VqS5A9fwmNmhwE/BP4Vbks8rScQB2aWBxxSw1O3O+emVGxzO/6ScGIi25YosXwGEWQ1PKYx2RFmZt8BXgauc859HXZ7QEUgLpxzPfb2vJldDvwcONul6cSMfX0GEbUW+H6179sCX4XUFgmZmTXGF4CJzrlXwm5PJXUHBczMegI3A+c75wrCbo8k1DzgKDM73MyaABcDr4fcJgmBmRnwFLDcOfdA2O2pTkUgeI8ALYDpZrbIzB4Lu0GJZma/NLO1wKnAW2b2TthtSoSKAQGDgXfwNwJfcM4tDbdViWdmzwMfAh3MbK2Z9Q+7TSE4HegLnFVxHFhkZueG3ShQbISISKTpSkBEJMJUBEREIkxFQEQkwlQEREQiTEVARCTCVARERCJMRUBEJMJUBEREIkxFQKSezOxgM/tfM/tjtcc6m1mhmf0qzLaJxEozhkUawMzOAd4AugGLgPnAXOfcb0NtmEiMVAREGsjMxgDnA+8DXYFc59yOcFslEhsVAZEGMrMsYDFwFHCacy4pFgsRiYXuCYg03GH4dQMc0D7cpojUja4ERBqgYqGQD4FV+OUC7wQ6O+dWh9kukVipCIg0gJmNBC4BOgPb8GtINwO6O+fKw2ybSCzUHSRST2bWDbgRuMw5t7Vi6dB+QCf8anIiSU9XAiIiEaYrARGRCFMREBGJMBUBEZEIUxEQEYkwFQERkQhTERARiTAVARGRCFMREBGJsP8HeMFB7Ah862wAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "\n",
    "# generate data\n",
    "np.random.seed(123)\n",
    "N = 20\n",
    "beta0 = -4\n",
    "beta1 = 2\n",
    "x = np.random.randn(N)\n",
    "e = np.random.randn(N)\n",
    "y = beta0 + beta1*x + e\n",
    "true_x = np.linspace(min(x), max(x), 100)\n",
    "true_y = beta0 + beta1*true_x\n",
    "\n",
    "# plot\n",
    "fig, ax = plt.subplots()\n",
    "sns.scatterplot(x, y, s = 40, label = 'Data')\n",
    "sns.lineplot(true_x, true_y, color = 'red', label = 'True Model')\n",
    "ax.set_xlabel('x', fontsize = 14)\n",
    "ax.set_title(fr\"$y = {beta0} + ${beta1}$x + \\epsilon$\", fontsize = 16)\n",
    "ax.set_ylabel('y', fontsize=14, rotation=0, labelpad=10)\n",
    "ax.legend(loc = 4)\n",
    "sns.despine()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Parameter Estimation\n",
    "\n",
    "The previous section covers the entire structure we assume our data follows in linear regression. The machine learning task is then to estimate the parameters in $\\bbeta$. These estimates are represented by $\\hat{\\beta}_0, \\dots, \\hat{\\beta}_D$ or $\\bbetahat$. The estimates give us *fitted values* for our target variable, represented by $\\hat{y}_n$. \n",
    "\n",
    "\n",
    "This task can be accomplished in two ways which, though slightly different conceptually, are identical mathematically. The first approach is through the lens of *minimizing loss*. A common practice in machine learning is to choose a loss function that defines how well a model with a given set of parameter estimates the observed data. The most common loss function for linear regression is squared error loss. This says the *loss* of our model is proportional to the sum of squared differences between the true $y_n$ values and the fitted values, $\\hat{y}_n$. We then *fit* the model by finding the estimates $\\bbetahat$ that minimize this loss function. This approach is covered in the subsection {doc}`s1/loss_minimization`.\n",
    "\n",
    "\n",
    "The second approach is through the lens of *maximizing likelihood*. Another common practice in machine learning is to model the target as a random variable whose distribution depends on one or more parameters, and then find the parameters that maximize its likelihood. Under this approach, we will represent the target with $Y_n$ since we are treating it as a random variable. The most common model for $Y_n$ in linear regression is a Normal random variable with mean $E(Y_n) = \\bbeta^\\top \\bx_n$. That is, we assume \n",
    "\n",
    "$$\n",
    "Y_n|\\bx_n \\sim \\mathcal{N}(\\bbeta^\\top \\bx_n, \\sigma^2),\n",
    "$$\n",
    "\n",
    "and we find the values of $\\bbetahat$ to maximize the likelihood. This approach is covered in subsection {doc}`s1/likelihood_maximization`. \n",
    "\n",
    "\n",
    "Once we've estimated $\\bbeta$, our model is *fit* and we can make predictions. The below graph is the same as the one above but includes our estimated line-of-best-fit, obtained by calculating $\\hat{\\beta}_0$ and $\\hat{\\beta}_1$. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 328,
   "metadata": {
    "tags": [
     "hide-input"
    ]
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYEAAAEeCAYAAABsaamyAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3dd3wU1RbA8d9NQiqhE3oXEKQbqiAEQkCQpiKCIqiogChFQIpPQURBQFABFX0qKoqKAj5LEkILTQggIlVAOgpJJJAQ0u/7YzYxhJRdsptJsuf7+eQDOzs7c3Z2d87MnTvnKq01QgghnJOL2QEIIYQwjyQBIYRwYpIEhBDCiUkSEEIIJyZJQAghnJgkASGEcGKSBIQQwolJEhBCCCfmlElAKTVcKaWVUrfl8PwMpVSRu4su0/tK/0tSSp1QSr2mlPI0O76CYNZnp5R6XCl1zLLNYwp6/UWZUirY8n191cHreUAp9a1S6rRS6rpS6qhS6nWllK8j11vYuZkdQCH1IRBsdhD5MBA4B/gCA4Cplv8/a2ZQBaTAPzulVFVgGbACeAxIKMj1F2VKqcFA8wJa3UTgDDAN4/fREpgBBCilOmit0woojkJFkkA2tNbnML4khZJSykNrnZjLLPu01sct/1+nlKoPPKGUGluQX3Qr4rQ7kz67+oArsFxrvTW/CzNju5lBKVUGWAiMB764hdefAj7RWs+w8iV9tNaRmR5vVkr9AywHugAbbI3BjrGZximbg/KSXZNC+jSlVH2l1I9KqTjLaeVLSimXLPM2V0p9r5S6bDnt3KaU6pRlntuUUp8ppU5a5vlTKfWuUqpsDuttopQKUUrFAV/b+Jb2Al5ABVvjtMw3WCl1RCmVoJT6XSnVVym1SSm1ydo4rdwmDZRSq5VSlyzrOqOU+kYp5WbN85njyLLcnkqpHZb1XlFKrVFKNcxhO+f5+WZ53SdA+nZYb1nGJ9au29bPVym1Tim1I5vpTZVSyUqpITm9thB6Aziotf6yIFaWJQGki7D8Wy19guW3mayUmpl5RsvvM1Yp5e+I+JRSdZVSnyilziujWfGSUmqjUqqEI9aXTpKA7VZjHDH0B9YAM4Fh6U8qpVoB24FywJPA/UA0EKaUujPTcqpiHLGOA3oArwDdgJ9yWO9aYDPQF+PoyRa1gSuWOGyKUynVHaOZ44hlnvnAIqCBtXHasE1+wPgxjsLYJlOARP79nub1/E2UUj2BH4E4YJDltU2ArUqpatm8JNfPNxuzgOcs/38GaG+ZZuu6rf18twKtlFIemd6jApYC27XWNx1RK4ObFX+uuazXrpRSHYFHgdEFtc4cdLb8ezh9guUs+kNgvFKqAoBS6iXgcWCA1nq3vYNQxlnRNqAMxvekG8Zv5RutdbK913cDrbXT/QHDAQ3clsPzM4xNc/M04LEs038HQjM9Xo/xhXLPNM3VMm1NLjG5AR0t62iZzXrH2vC+GlqWVxbji5sCjMkyr1VxYuy8DwAq07RWlvVssiZOa9aFcZaigb45vLdcn8/pswN2A8cAt0zT6gDJwJu2fr45rDPQ8touWabnuW5bPt8s62qXadowyzKb5vCaLpbX5PW3yZoY7PD7KwEcBF7NNE1nfpzNa5TlO5357xTGwVPmaa42xFENuASsy+a5ysA1YB7wBJAKPOio2IBelm0wGPDO/J1x+OdRUCsqTH/kLwn4ZZn+JXDE8n8vjB1u1g/fDXgH+CfT69wxLlAdAa5n+TE+lM16a9rwvrL+Lckyn1VxYuyok4CZ2azrT7JPAjVvcV0KOAEcwjgCqp9lObk+n91nB/gAaWSzc8Fowtljy+eby3a/KQlYu25bPl/L/CUt23Oc5XEZ4CKwMJfX+AL+Vvw1tGL96e/1lhMK8KLl++OVaVpeSaBLftebzXbcDVwAqucwz2yMi/wpwDOOjA2oCPyRad5oa96HPf7kwrDt/snyOBFI735ZDmPH+R/L302UUi7auDj7OkZvnVcwjrZjgerAd5mWl9lfNsQ4AKOpqSIwARitlNqptf7Uljgxjr5LYBwtZXUxh3VnjdPqbWJpepqBsW3KK6VOAvO01u9qrXVuz+cQS1mM5JHdtvsbqJXN9Nw+X1vYum6rPl+tdZxS6jegnWXSbIxk83IuL4sD9lmzeCvm2Q40smK++OwmKqVqAtOBEYBH5mYty+MyQKzWOjXLS/cArbNM+x6jiXBZpmmxeQWmjO7S3wN1gc7a6EyQnWOAB7BVa70kl0XaIzYf4GfgTWA/UGDdjCUJ2FcMxg9yCfBpdjPof3vnPAR8qrXO6ButlCqZy7Kt+YGmO6AtvYOUUhswvlTzlFLfaq2vWRunUioKo5nBL5tZKmF0t8srTqu3idb6T+BRSxt3c2AMsFQpdUpr/XNez2ez6MuWeCpn81xlMl0jcQBb123L57sN6Gu51jISGKa1vprL/J2BjVYsdzPGUW2OtNbxGGevt6ouRlL9PJvnJlr+WpIlaWmtYzGO3DMopZKAC9qGNnrLRdZvgTZAoNb69xzm6wq8D+wA7lJKNdda/5bdvPmNTSlVCuPs8Fmt9f+sfS/2IknAjrTW15RSWzB2UHt17t0xvTF2sJk95oCYEpVSkzAuPI7GOHK2Kk6tdapSajdwv1Jqhract1ou5tYh+ySQdRm2bJP012hgn1JqAkZ7bBOMoySrns+y7j3AQEv8qZb4awEdMJqjHMLB696KcRb5KbBNa53dDjWz7I5Us5PnUbQd7AMCspm+ESMx/Bc4ns3z+WY5u12BcdG1t9b6lxzma4XRKeBDjO6rfwCvAb0dERdwN8aZ4VEHLT9Xzp4Eeiql/s4y7Uo+lzkBCAdClFL/xTjNr4BxMdVVaz3FMl8wMEwp9TvGl/4+jJ2D3Wmtv1dKRQATlVKLtdbXbYjzZSAUWK2UWmaZZwZGk4a19xzkuS6lVDPgLeArjO3hinGNIwXYkNfzuaz7Pxg9dH5QSi3FaAueifE5L7Ay/lvlqHVvs/x7O8Y2zFV2R6pm0VrH8G+X2gzGyR2ntdY3PWdHSzBupJwNXFNKtcv03Dmt9TllVBH4GeM7/6zljHgm8JFS6m6tdbgD4krvuvqpUmoexllidaCr1vpxB6zvRgV18aEw/ZHzBVSN0RNmBjlfGHbLMv0T4FSWaY2AlRht6YkY7fPfA70yzVPBMs9ly98KjKM1DQzPa715vK+bLngDQZbnxtsSp2W+IRhHKYkYvToGAL8Cq62NM691YTQ5Lcc46orHaJvfDPSw5vmscWSZ1hPjtP46xg54LVkugtry+Wbz3rLtHWTNum35fDO9prRlGy4y+7dkx99krheGc3jNKWCGjfPn9LufgdFM9ydGkvLI9Lr0nmzbHRjbIxhnbFct3+9DwIsFse2VJQAhrKaUqo5xND5baz3L7HicjVJqAUZivl1rnd8zV+HknL05SORBKeWF0WMhDIjCuLA3GeNo5UMTQ3MqSilvjOsqnYCxwEBJAMIeJAmIvKRinCYvBspj3ECzBWMnZEu3VZE/gRhNSecxbixbbXI8opiQ5iAhhHBiUjtICCGcWJFqDurZs6cODi7KZf6FEMIUKqcnitSZQFRUlNkhCCFEsVKkkoAQQgj7kiQghBBOTJKAEEI4MUkCQgjhxCQJCCGEE5MkIIQQTqxI3ScghBDOIj4xhdjEFC7HJ1HW2x1fDze8Pey/y5YkIIQQhUx8Ygphhy8y8Zv9JKWm4e7qwvyBzQhsVMnuiUCag4QQopCJTUzJSAAASalpTFq1n9jEFLuvS5KAEEIUMpfjkzISQLrElDRi4rOOSJt/kgSEEKKQKevtjrvrjbtnDzcXyniXsPu6JAkIIUQh4+vhxvyBzfBwM3bRHm4uzHugOb4OuDBcpMYT8Pf317t3F4rxsoUQwqHSewfFxCdTxrtEfnsH5VhFVHoHCSFEIeRt2elXKuXp0PVIc5AQQjgxSQJCCOHEJAkIIYQTkyQghBBOTC4MCyFEIVUQ9YMkCQghRCFUUPWDpDlICCEKofT6QX7RF7jj4gmH1Q+SJCCEEIVQzOVYntryBWH/Hc3skCWgtUPqB0kSEEKIwiYsjHqBdzFxy+eE1WvDyP7TQCmH1A+SawJCCFFYXLgAEybAV1/hUu82diz5nOcvlCMxJc1h9YMkCQghhNlSUuCdd+DllyEpCWbOxGXyZJorN8LtVz8oW5IEhBDCTNu2wejRsH8/3HOPkQzq1QPAGxxeP0iuCQghhBkiI+GJJ6BjR7h8Gb77Dn78MSMBFBRTk4BSqqdS6qhS6rhSaoqZsQghRIFIS4Nly6BhQ/j0U3jhBTh8GAYMAJVjxWeHMS0JKKVcgSXAPUBjYLBSqrFZ8QghhMPt3QsdOsDTT0OzZvDbbzBnDvj4mBaSmWcCbYDjWus/tdZJwEqgn4nxCCGEY1y5As89B61bw8mT8NlnsHEjNDb/uNfMJFANOJvp8TnLtBsopZ5SSu1WSu2OjIwssOCEECLftIYVK4ymnyVLYNQoOHoUHnnE6qYfrTVJcUkOC9HMJJDdFrhprEut9TKttb/W2r9ixYoFEJYQQtjB4cPQrZuxw69ZE3btgsWLoUwZqxdxPuI8H931Ed8/8b3DwjQzCZwDamR6XB24YFIsQghhH9euwdSp0Lw57NsH770HO3bAnXdavYjYC7GsGbaGD9t8yOU/L1Ovp+N6DJl5n0AEUF8pVQc4DzwEDDExHiGEuHVaw/ffG23/Z87A8OEwdy74+Vm9iOTryex4cwdbX99KWnIad025i07TOuHh6+GwsE1LAlrrFKXUGCAEcAU+0lofNCseIYS4ZSdPGjv/H36AJk1gyxaj/7+VtNYcWnWIdZPWceX0FRrd14jANwIpV6+cA4M2mHrHsNb6J+AnM2MQQohblpgI8+bB7Nng5gYLFsCzz0IJ64u8/fXrXwSPDebMljNUalaJfhv6USegjgODvpGUjRBCiFuxbh2MGQN//AEDB8Kbb0L16la/PO7vONZPX8++j/fhXcGbe9+/l5ZPtMTFtWAv1UoSEEI4VEEMkVigzp+H55+Hr76C226D4GDo0cPql6ckprDzrZ2EvxpOyvUU2k9oz93/uRvP0o6rD5SbIvxJCCEKu4IaIrFApFf6fOklSE6GmTNh8mTwtG7nrbXm6NqjhE4M5fKJyzS4twFBC4Io36C8gwPPXRH7FIQQRUn6EIlJqWkAGUMkhk8OKFpJYNs240av33+HXr3g7bdtKvR2cf9FQsaHcHLDSSo2rsgjIY9QL6hgC8XlpAh9CkKIouZyfFJGAkiXPkSiI8sj201kpFHg7eOPoUYNo9Jn//5W3+17LfIaG/+zkb0f7MWzjCf3LL4H/6f9cXErPAWcJQkIIRymrLc77q4uNyQCRwyRaHdpafDBB8ZNX7GxRrPPSy9ZXegtNSmVXUt2sXnmZpLikmg9pjVdXu6CVzkvBwduO0kCQgiH8fVwY/7AZkxatd+hQyTa1d69RtPPrl3QpYtR88fKQm9aa479eIzQ50OJ/iOa23reRtCbQVRsVHhL3hTiT0IIUdR5e7gR2KgS4ZMDHDpEol3ExMB//gNLl0LFikalz4cftrrpJ/JQJCETQjgRcoLyDcsz5Mch1O9V38FB518h/CSEEMWJt2WnX2ivAWgNX3xhdPuMjDSGepw1y+pCb/HR8WyasYnd7+7GvaQ7PRb2oPUzrXEt4ergwO1DkoAQwnkdPmzs9DdtgjZtjOEdrSz0lpqcyu73drPp5U0kXknkzqfvJOCVALwreDs2ZjuTJCCEcD7XrsGrr8L8+eDra1T6fPJJcLGu187xkOOEjA8h6nAUdbrVocfCHlRqWsnBQTuGJAEhhPPQGtauhbFjjUqfjz1mVPq0cqySqKNRhD4fyrEfj1G2XlkGrRlEw74NUSaMDWwvkgSEEM7hzz+NSp8//mhzpc+EmAQ2v7KZXe/sws3LjcA3Amn7XFvcCuMFbhsV/XcghBC5SUyEN96A116zudJnWmoaez/Yy8b/bCQ+Op6WT7Sk66tdKVmpZAEEXjAkCQghiq916+CZZ+DYMZsrfZ7ccJLgccFc+v0Ste6uRY9FPajSsoqDAy54kgSEEMXP+fMwYQJ8/bXNlT7/OfEP6yau48iaI5SpXYaB3wyk0f2NinS7f24kCQghio/kZKPS58svG1U/baj0mXg1kfDZ4exctBOXEi50nd2V9hPa4+ZZvHeTxfvdCSGcx9atRp//9Eqf77wDdevm+bK01DT2fbKPDdM3cO3iNRoMacLtk9pTtW45klTx30kW9/cnhCjuslb6XL0a+vWzqtzD6fDTBI8L5u9f/6Zqu+pUWxDIy4f+Imnlr0V77AMbFJ56pkIIYYu0NHj/fWjY0KjzM3mycQewFaWeY07F8M3Ab/ik8yfER8Vz/5f30yd4iJEAsox9EJuYUhDvxjTFN70JIYqvPXuMph8bK30mxSWx5fUt7FiwAxdXF7rM7EKHiR0o4V2CI39fLdpjH9wiSQJCiKIja6XPzz+HIUPyPPLXaZr9n+8nbEoYcX/F0fThpgTOCaRU9VIZ8xTZsQ/ySZKAEKLw0xpWrICJE22u9Hl2+1mCxwVzIeIC1dpUY9B3g6je7uZ7BYrk2Ad2ULzfnRCi6Dt0iNSRo3DdEs71lncS//VqvNq2zvNi7ZWzV1g/ZT2/f/E7vlV96f9pf5o93Azlkv1ZQ5Ea+8COive7E0IUXdeuwaxZ6AULSPX2YcY9z/J50+6UCIlhfqmLOfbaSY5PZtu8bWybuw00dHqxEx1f6Ih7Sfc8V1noxz5wAEkCQojCJb3S53PPwdmzJDzyKF3K38NFT1/g31474ZMDbkgCWmsOrDxA2AthXD17lTsevIPANwIpU8u6wWGclSQBIUThkbnSZ9Om8OWXnK7XlIuLttwwW9ZeO+cjzhMyLoSz289SpVUV7ltxH7U61TLjHRQ5kgSEEObLpdJn2asJOfbaib0Qy/pp6/lt+W/4VPKh73/70nxYc1xc5RYoa0kSEEKYKzQUxowxKn0++KBR6bNatYyns+u1M7fPHRxctJNf3thGWnIad71wF52mdcKjlIeJb6RokiQghDBH5kqf9etDSAgEBd00W+ZeO5evJfFP6J/sfPA7rp65wu39b6f7/O6Uq1fOhDdQPJiSBJRS84A+QBJwAnhMax1jRixCiAKWtdLnrFkwaRJ45HwU7+3hxpVDkewYG8yZLWeo1KwS/T/uR52udQow8OLJrDOBdcBUrXWKUmouMBV4waRYhBAFJXOlz9694e2386z0GXcxjg3TN/DrR7/iXd6b3u/1ptWIVtLubyemJAGtdWimh78AD5gRhxCigERGGgXePvkEataENWugb99cyz2kJKaw862dhL8aTsr1FNqNb0fn/3TGs4zz9OEvCIXhmsDjwFdmByGEcIDUVPjwQ5g6FWJjYcoUePFF8PHJ8SVaa46uPUroxFAun7hMgz4NCJofRPkG5QswcOfhsCSglAoDKmfz1HSt9VrLPNOBFGBFLst5CngKoGbNmg6IVAjhELdQ6fPi7xcJGRfCyQ0nqdi4Io+EPEK9oHoFE6+TUlprc1as1DBgJNBNax1vzWv8/f317t27HRuYECJ/YmKMo/133zUqfb75JgwenGvTz7XIa2x8aSN7l+3Fs4wnXWZ2wX+kPy5u0u5vJzlufLN6B/XEuBDc2doEIIQo5NIrfT7/PERFWVXpMzUplV1LdrF55maS4pJo/Uxruszoglc5rwIM3LmZdU1gMeABrFPG0cEvWuuRJsUihMivQ4fgmWdg0yZo0wZ+/hlatcpxdq01x346RujzoUQfjaZej3r0eLMHFRtXLLiYBWBe76DbzFivEMLOLJU+WbAAfH2N4R5HjACXnJtxIg9FEjIhhBMhJyjfoDyDfxhM/V71UVaMCSzsrzD0DhJCFDVZKn3y+OMwZ45xDSAH1/+5zqYZm4hYGoF7SXeC3gyizTNtcHV3LcDARVaSBIQQtsmm0id33ZXj7KnJqex5fw+bXt5EQkwCrZ5qRcArAfhUzLmbqCg4kgSEENbJWunzzTeNSp9uOe9GjoccJ2R8CFGHo6jTtQ49FvagUrNKBRi0yIskASFE3jJX+hw0yLgGkKnSZ1bRf0QT+nwof/zwB2XrlWXQ6kE07NdQ2v0LIUkCQoicnT8P48fDN9/kWukzXUJMAptf2cyud3bh5uVG4NxA2o5ti1sxH6e3KJNPRghxs6yVPl95xaj9k0Olz7TUNPZ+uJeNL24kPjqelk+0pOurXSlZqWQBBy5sJUlACHGjzJU+e/UykkEulT5PbjhJ8LhgLv1+iVp316LHoh5UaVmlAAMW+SFJQAhhsLHS5z8n/mHdpHUcWX2EMrXLMPCbgTS6v1Gha/ePT0whNjGFy/FJlPV2x9fD7YYB6p2dbAkhnF1qKnzwAUybZlWlz8TYRLbM3sIvC3/BpYQLXWd3pf2E9rh5Fr7dSXxiCmGHLzLxm/0kpabh7urC/IHNCGxUSRKBhWwFIZzZnj0wahRERORZ6TMtNY19n+xjw/QNXLt4jebDmtPttW74VvXNdRVmHonHJqZkJACApNQ0Jq3aT/jkAEkCFrIVhCjibmknm17pc+lS8PODzz6Dhx/Osenn9JbTBI8N5u9f/6Z6++oM/n4w1drk3EU0c2xmHolfjk/KSADpElPSiIlPplIpGZwGJAkIUaTZvJPNWulzzBij508OlT5jTsWwbvI6Dn1ziFLVS3HfF/fR5KEmVrf7m30kXtbbHXdXlxsSgYebC2W8Szh83UWFFOsWogjLaScbm5hy88wHD0JAAAwdCrVrG01Ab7+dbQJIiktiw4sbWHz7Yv744Q86z+jMmKNjaDq4qU0XfnM7Ei8Ivh5uzB/YDA/LuAQebi7Me6A5vtIUlEG2hBBFmFXNHXFxRqXPN980Kn0uWwZPPJFtpU+dptn/+X7WT11P7IVYmg5pSrc53Shdo/QtxWf2kbi3hxuBjSoRPjmAmPhkyniXkN5BWciWEKIIy3Unq7XRzXPsWKsqfZ7dcZaQcSGc33Weqq2rMnDVQGq0r5Gv+NKPxCet2k9iSpopR+Lelp2+XAPInmnDS94KGV5SiBulXxPIupPt7h6L18TxxuAuTZsaQz3mUOnzytkrrJ+ynt+/+J2SVUoSOCeQZo80Q7nYp79/+oVrORI3VY4fpiQBIYq4G3ayLqmUXbwI9zfmQIkSxkXfHCp9Jscns23eNrbN3QYa2k9sT8cXOuJe0t2EdyEcrHCNMSyEsJ+M5o5fwo0hHo8fhwcfNK4BZFPpU2vNgZUHCHshjKtnr9J4YGO6v9GdMrVzHgtYFF+SBIQo6s6dgwkT/q30GRoK3btnO+uF3RcIHhvM2e1nqdyyMvd9fh+17q5VwAGLwkSSgBBFVdZKn7NmwaRJ2Vb6jL0Qy/pp6/lt+W/4+PnQ58M+tBjeAhdX6SXu7CQJCFEUZa702bu3kQzq1LlptpSEFHa8uYMtr20hLTmNDpM7cPf0u/EolX1JaOF8JAkIUZRYWelTa83hbw+zbtI6Yk7FcHv/2+k+vzvl6pUzJ25RaOWZBJRSjwILgapa68RM01cAvlrrvg6MTwgBRqXPDz+EqVONm79yqfT5169/ETIuhNPhp6nUrBKPrn+UOl1vPksQAqw7E/gGeAvoB3wNoJQqDQwABjsuNCEEYFT6HD0adu0yKn0uXQqNGt00W9zFODa8uIFf//sr3uW96f1eb1qNaCXt/iJXeSYBrfV1y1H/41iSADAEuAr86MDYhHBu6ZU+333XuMt3xQoYPPimpp+UxBR2vrWT8FfDSbmeQrvx7ej8n854lpE7ZEXerL0m8AGwVylVXWt9DiMhLNdaZ1OlSgiRL1rD55/DxIlGpc/Ro42eP1kKvWmtObr2KKETQ7l84jIN7m1A0IIgyjcob1LgoiiyKglorX9TSu0Fhiul1gD+wCMOjUwIZ3TwoLHTDw+Htm2Nsg+tWt0028XfLxIyPoST609SsXFFHgl5hHpB9UwIWBR1tvQO+gCYDFQAtmmtjzomJCGcUOZKn6VK5VjpMz4qno0vbWTP+3vwKO3BPe/cg/9If1zcpN1f3BpbksCXwJvAKGCkY8IRwsloDatXG5U+z50zdvxz5kCFCjfMlpqcSsSSCDbP3ExibCKtn2lNlxld8CrnZVLgoriwOglorWOVUl8DA/n3ArEQ4ladOAHPPQc//QTNmsHKldlW+jz20zFCJoQQfTSaekH16LGwBxUbZ18OWghb2XqzWBVgpdb6miOCEcIpJCTAG2/Aa68ZlT4XLjSGecxS6TPycCShE0I5Hnyc8g3KM/iHwdTvVd+mkb2EyItVSUApVQ4IBIKA5vYMQCk1EZgHVNRaR9lz2UJkdksDsttbSIixwz9+HAYNggULbqr0ef2f62yauYmIJRG4l3QnaEEQbca0wdXdtWBjFU7B2l/AXqAcME1rfcBeK1dK1QC6A2fstUwhsmPzgOz2du4cjB8Pq1blWOkzLSWN3e/tZtPLm0iISaDVk60ImBWAT8Wb7woWwl6s7SJa20HrX4jR42itg5YvBJDzgOzhkwMcmwSSk43B3GfMMCp9vvqq0f8/S6XPE6EnCBkfQuShSOp0rUOPhT2o1KyS4+ISwsK0AnJKqb7Aecs9CLnN9xTwFEDNmjULKDpR3Fg1ILu9bd0Ko0bBgQM5VvqM/iOa0OdD+eOHPyhbtyyDVg+iYb+G0u4vCoxDk4BSKgyonM1T04FpGNcYcqW1XgYsA2N4SbsGKJxGrgOy29ulS/DCC7lW+kyISWDzrM3semcXbp5uBM4NpO3YtrjJ2LuigDn0G6e1DsxuulKqKVAHSD8LqI5RlqKN1vpvR8YknJOvhxvzBza7aUB2X3vudFNT4YMPjEqf165lW+kzLTWNvR/uZeOLG4mPjqfl4y3pOrsrJSuVtF8cQtigUAw0r5Q6Bfjn1TtIBpoX+XHDgOzeJezbO2jPHqPpJyICAgJgyZKbKn2e3HiSkHEhXNx/kZqdatLzrZ5UaVnFPusXIncy0LwQGQOy2/MaQHqlz6VLwc/PKPw2ZMgNTT+X/7xM6MRQjqw+QulapXng6wdo/EBjafcXhWZauuEAAB8sSURBVEKhSAIO7H0khGNkrfQ5ZoxR+6d06YxZEmMT2TJ7C78s/AWXEi4EvBpA+wntKeFl+3WIQnGPgyiW5FskhK3yqPSp0zT7PtnH+mnruXbxGs0fbU6317vhW9X3llZn+j0OoliT0oNCWCsuzuj106KF0e1z2TLYvv2GBHBm6xk+aP0B3z/xPWXrlmXEzhH0X97/lhMA5HyPQ2yiDOch8k8OI4TIixWVPmNOxxD2QhgHvzpIqeqluG/FfTQZ3MQu7f6m3OMgnIYkASFyc/y4Uenz55+NSp9ffQUdOmQ8nRSXxNa5W9kxfwco6PxyZzpM6oC7j7vdQijQexyE05EkIER2EhJg7lx4/fVsK33qNM3+FftZP2U9sRdiaTqkKd3mdKN0jdJ5LNh2BXKPg3BaheI+AWvJfQKiQAQHGzv8EyeyrfR57pdzBI8N5vyu81RtXZWeb/WkRvsaDg3Jofc4CGcg9wkIkafMlT4bNIB16yDw35ver567StgLYfz+xe+UrFKS/sv70+yRZigXx/f3d8g9DkIgSUAIo9LnW28ZlT5TU43+/pMmZVT6TI5PZvv87WydsxWdpuk0vRMdp3TEvaT92v2FMIskAeHctmwxyj0cPAj33muUfbZU+tRac/Crg6ybvI6rZ6/SeGBjur/RnTK1y5gctBD2I0lAOKdLl2DyZFi+HGrVgrVrjUqfFhd2XyB4XDBnt52lcsvK3Pf5fdS6u5aJAQvhGJIEhHPJWulz6lSYPj2j0mfsX7FsmLaBfZ/sw8fPhz4f9KHFYy1wcZX7KkXxJElAOI9cKn2mJKSwY+EOtr62ldSkVDpM7sDd0+/Go5RHHgsVomiTJCCKv5gYmD4d/e67pFX04+LiD3B95GF8PUvgpTWHvz3MuknriDkVw+39b6f7vO6Uu62c2VELUSAkCYjiK1OlTx0VxclBw3mgai/+OeuF+xubeKVpVa4t2cu5rWfwa+rHo+sfpU7XOnkvV4hiRJKAKJ6yVPr8Z9VaegZfJik1Dc9rybQKP8+R13bhVd6L3u/2ptWIVri4Sbu/cD7yrRfFS1yc0esnS6XPyPqNSUlMocnOv7l/2e/UPxDNQf9K9Nr2GP4j/SUBCKclZwKieNAavvsOxo27qdKn1pqYjacZ8NFBSl1O5Ey90kR0rUGinzd+VW69xLMQxYEkAVH0nThh1PoJDobmzW+o9Hnx94uEjA/h5PqT+NUpw889anGqVikpwiaEhfwCRNGVudKnuzssWgTPPANubsRHxbPxpY3seX8PHqU9uOede2g0vDnD07QUYRMiE/kFiKIpc6XPhx4yKn1WrUpqcioRb/3C5hmbSYxNxH+0P11mdMG7vDcAviBF2ITIRJKAyLcCHQT97Fmj0ue3395U6fPYT8cImRBC9NFo6gXVI+jNIPzu8HNMHEIUE5IERL4U2CDoWSt9zp4Nzz8PHh5EHYkiZHwIx4OPU65+OQb/bzD1e9e3y9COQhR3kgREvuQ0CHr45AD7JYEcKn1e/+c6myb/TMSSCNxLuhO0IIg2Y9rg6u5qn/UK4QQkCYh8cegg6FkqfSau+o6Y7vcQfTWBv+ZvY/fr20iMSaDVk60ImBWAT0Wf/K1PCCckSUDki0MGQU9NNW7ymjYto9Jn/MQXCDsTx/yn19Bq3WnKRidQsk1VBi3pRW3/ankvUwiRLblNUuRL+iDoHpY7bvPd/37PHmjf3ij50KoV7N8Pr73G2VPX+N+Dq+i28iiuqZr1A+qxLLAaXg3K2/HdCOF85ExA5Iu3hxuBjSoRPjkgf/3vL1+GF1+Ed9+FSpVgxQoYPJiEK4lsfj6EXe/swk9BRJdqHLqzEmluLpCq7dPsJIQTkyQg8i1fg6BnqvRJVBQ8+yy88gppJX3Zu2wPG1/cSHx0PLc/0ozZZV246vXvVzbfzU5CCEkCxVGB9tvPjyyVPgkOhpYtObXpFMHjvuTibxep2akmPRf1pPQdFfE6fJFJq/aTmJImZR+EsBP5BRUzBdZvPz/i4uCVV2DhQihVyhju8fHHuXzqCqH3fcWR1UcoXas0D3z9AI0faJzR398uzU5CiBvIL6iYKZB++7cqh0qfiR6+bJm+gV/e/AWXEi4EvBpA+wntKeF1Y1NPvpqdhBDZMm2voJR6FhgDpAA/aq0nmxVLceLQfvv5cfy40d6fqdKnbteefcv3sWHaBuL+jqPZ0GZ0e70bpaqVynNxRabJS4hCzpRfjVIqAOgHNNNaJyqlpMCLnTik335+ZK30+dZbMHo0Z365QHCbD/hrz19Ub1+dh9Y+RLU21vX3LxJNXkIUEWbdJzAKmKO1TgTQWl8yKY5ix+799vMjOBiaNDHq/QwYAEeOENPvUVY9vIaPO33MtYvXGPD5AB7f9rjVCQBybvKKTUxx0BsRovgy67CpAdBJKTUbSAAmaq0jsptRKfUU8BRAzZo1Cy7CIspu/fbzI3Olz4YNISyMpLad2Dp3Kzvm7wAFnV/uTIdJHXD3cbd58YW2yUuIIshhewalVBhQOZunplvWWxZoB7QGvlZK1dVa66wza62XAcsA/P39b3pe3My0C6jJycbALjNnQloavPoqesLz7F91lPUNFxN7IZamQ5rSbU43StcofcurKXRNXkIUYQ5LAlrrwJyeU0qNAr6z7PR3KaXSgApApKPiEQ4WHm70+T94EPr0gbff5tzfbgR3+Zzzu85T1b8qA78ZSI0ONfK9qvQmL7lnQIj8M+tXswboCmxSSjUA3IEok2IR+XHpEkyaBJ9+CrVqwdq1XG3VhbApYfy+4ndKVilJv0/60Xxoc5SLfer756fJS3oVCXEjs779HwEfKaUOAEnAsOyagkQhlpoK778P06cblT6nTSN5/GS2L93HtsGLSUtNo+O0jnSa2gn3kra3++flVpq8pFeREDcz5ZuvtU4CHjFj3cIOIiKMQV727IGuXdGLF3Pwt1TC7vyYK2eu0HhgYwLnBlK2TtkCC8maI/xCfSOdECaRb76w3uXLxpH/e+8ZlT6//JIL9ToR/GQIZ7edpXKLyvT/tD+1O9cu0LCsPcKXXkVC3EzGExB509po82/Y0GgCeu45YjftYW2IFx+0+ZB/jv1Dnw/68OTuJws8AYD19w2k9yrKTHoVCWcnZwIidwcOGL1+tmyBdu1I+d/P/LLxOlv8l5OSmEKHSR3oNL0TnqXNO5K29ghfehUJcTP59ovsZan0qZct43CZDqwbHEbMyRga9mtI0Pwgyt1W7pZXYa+eOtbeN1AobqQTopCRb7+4UTaVPv9++HmCZ+7k9OZV+DXxY2jYUOp2q5uv1dizp44tR/hSiVSIG6mi1DPT399f79692+wwiq8slT6vzV7EhrWx7P1wL97lvQmYFUCrEa1wccv/paSLVxPoNHfjTUfv4ZMDbmkHnX5WIUf4QmQrx5t05Fcibqr0mTp/ITtT/QkfspXk+GTajWtH55c641nGfkfP9u6pI0f4QtwaSQLOLjgYxoyBEyfQDw3mjx7PEjp7N/8cX0/93vUJWhBEhYYV7L5aqf8jROEgScBZZan0eemDtYR8dYU/HwulQqMKPBz8MLf1uM1hq5eeOuZJTk7m3LlzJCQkmB2KsDNPT0+qV69OiRLWH0zJNQFnk6XSZ/yE6WyMasaeD37Fo7QHXWZ2wX+kP64lXB0eirTjm+PkyZP4+vpSvnz5jPGbRdGntSY6OprY2Fjq1KmT9Wm5JiC4odJnau8+RNw5ks1v/05i7K/4j/any4wueJf3LrBwpB3fHAkJCdSuXVsSQDGjlKJ8+fJERtpWjFmSgDO4eBEmT86o9HnspU8J/foqUT9GULd7XXos7IHfHTLC560qipVJJQEUT7fyuRbub6rIn9RUWLYMpk2Da9eIemoaIacacvyVPylXvxyD/zeY+r3ryw4hH6QyqSjqpHZQcRURAe3awejRXG/WluCHPuHdjzw5u/MCQQuCGH1gNA3ubSAJIJ9kvGPbRUdH06JFC1q0aEHlypWpVq1axuOkpCS7rScsLAylFMuXL8+YFhERgVKKRYsWWb2c48eP06JFi3zPU1jJoUpxk6nSZ5pfZXYPX8qm72NJ2HqcliNa0nVWV3z8fMyOstiQyqS2K1++PPv27QNgxowZlCxZkokTJ94wj9YarTUuLvk7Tm3atCkrV65k2LBhAKxcuZLmzZvna5nFjSSB4kJr+OwzmDgRoqM50W8CIUdrE/nJJWoH1Kbnop5UalbJ7CiLnSJ/v8O4cWDZIdtNixZGDzQbHT9+nP79+9OxY0d27tzJmjVraN68OTExMYCxAw8LC+PDDz/k4sWLjBo1ijNnzuDi4sLbb79Nu3btblpm3bp1iYyMJCoqivLly7Nu3TruueeejOf37t3LqFGjuH79OvXr1+ejjz6idOnSRERE8MQTT+Dj48Ndd92VMX9KSgqTJ09m69atJCQk8NxzzzFixIhb2EiFhzQHFQcHD0LnzjBsGNFVm/Jlx6V8vsaXlMQ0HvzuQR5d/6gkAAdJv9/Bw1JKQ+53yJ9Dhw7xxBNP8Ouvv1KtWrUc53vuueeYPHkyu3fv5uuvv851R3z//fezatUqwsPDadu27Q196B955BEWLFjA/v37adiwIbNmzQJg+PDhvPvuu+zYsYPU1NSM+ZctW4afnx+7du0iIiKCJUuWcObMGTu8c/PIN7Uoi4sz+vsvWkRCyQqEB81n58Z43Dz+IXBuIG3HtsVNdkYOVeQrk97CEbsj1atXj9atW+c5X1hYGEePHs14fPnyZa5fv46Xl9dN8w4aNIihQ4fSoEEDBg8ezIYNGwDj2kRCQgIdO3YEYNiwYQwdOpSoqCiuX7+ecQYwdOhQNm7cCEBoaCiHDx9m5cqVAFy5coVjx45Rq1at/L1xExWRb6q4QaZKn2nnzvNrx+fYcKQq8eviaPFYC7rN7kbJyiXNjtJpyP0O9uPj8+/1KhcXFzLfzJr5DmetNbt27cLdPe/xq6tVq4bWms2bN7N06dKMJJDbjbI5dZjQWrN06VK6det2w/Tjx4/nGUdhJc1BRc3x49CrFzzwAKc8GrKs/nx+2FqWCrdX5KndT9Hvv/0kAYhiwcXFhbJly3Ls2DHS0tJYvXp1xnOBgYEsWbIk4/G+PK5rzJo1i7lz595woblChQp4eXmxfft2AD777DM6d+5MhQoV8PT0ZMeOHQCsWLEi4zU9evRg6dKlpKQYvb+OHj3K9evX8/9mTSRnAkVFQgLMmQNz5nDZrSLrms3m8P5kStd04YGvHqDxwMbS3VMUO3PnzqVnz57UrFmTxo0bk5iYCMCSJUsYNWoUH3/8MSkpKQQEBNyQFLJKb/LJ6rPPPsu4MHzbbbfx8ccfA/Dxxx8zYsQIfHx8CAoKypj/6aef5syZMxndQf38/Fi7dq293q4ppHZQUfDzz/DssySeOMuWRk/zy4mKuLi50nFqR9o/354SXkWkJ4ooFA4fPkyjRo3MDkM4SA6fr9QOKpLOnoVx49DfrWZf5R5sKPc4cYeTaTa0Cd1e70apaqXMjlAIUcRJEiiMMlX6PJNcheCqM/jrgqZ6u0oMWtSD6m2rmx2hEKKYkCRQ2FgqfcYcPEdY1Sc4eKEcvqokAz4PpOmQptLuL4SwK0kChcXFizBpEkmfrWRb6d5sd38Q/nHl7pc6cNfku3D3ybsrnBBC2EqSgNlSU+H999FTp7H/Wl3Wl5xG7BVFk8F3EDgnkNI1S5sdoRCiGJMkYKaICBg1inN7LhJc6gnOp5ai6u1VGfhWT2p0qGF2dEIIJyA3i5nh8mUYPZqrbQL57uDt/JcRXPGpSr+P+zFi5whJAKLYc3V1zSgf3aJFC+bMmZPjvGvWrOHQoUMZj1966SXCwsLyHUNMTAxLly61+XUzZsxg/vz52U5XSt1w9/DChQtRSmFL1/ZPPvmEMWPG5Hsea8mZQEHSGj79lOSJU9ge3ZBtbuNJ0250nNaeTlM74V5S2v2Fc/Dy8srzLt90a9as4d5776Vx48YAvPLKK3aJIT0JjB492i7Lg39LV7/44osArFq1KiPuwsopkkChGP7vwAH0qNEc3PoP69yHcVV70bh/YwLfCKRsnbIFG4sQFsHjgvl73992XWblFpXpuajnLb12ypQpfP/997i5uREUFMR9993H999/z+bNm3n11Vf59ttvmTVrFvfeey8PPPAAtWvXZsiQIWzcuJHk5GSWLVvG1KlTOX78OJMmTWLkyJHExcXRr18/Ll++THJyMq+++ir9+vVjypQpnDhxghYtWtC9e3fmzZvHvHnz+Prrr0lMTGTAgAHMnDkTgNmzZ/Ppp59So0YNKlasyJ133plt/P3792ft2rW8+OKL/Pnnn5QuXfqGqqVffvklr732Glprevfuzdy5cwHjDuXXX3+dKlWq0KBBAzw8PACIjIxk5MiRGZVKFy1adENpa3swLQkopVoA7wGeQAowWmu9y97rMX34P0ulzwtvfkmwS2/OUpXKjSsxYFFPaneu7fj1C1EIXb9+/YaRuKZOnUr37t1ZvXo1R44cQSlFTEwMZcqUoW/fvhk7/ezUqFGDHTt2MH78eIYPH862bdtISEjgjjvuYOTIkXh6erJ69WpKlSpFVFQU7dq1o2/fvsyZM4cDBw5knJGEhoZy7Ngxdu3ahdaavn37Eh4ejo+PDytXruTXX38lJSWFVq1a5ZgESpUqRY0aNThw4ABr165l0KBBGaUoLly4wAsvvMCePXsoW7YsQUFBrFmzhrZt2/Lyyy+zZ88eSpcuTUBAAC1btgRg7NixjB8/no4dO3LmzBl69OjB4cOH7flRmHom8AYwU2v9s1Kql+VxF3uvJKfh/8InBzg2CVgqfcaNmcL6v+9gHyPwKe9Nn9cCafFYC1xc5XKMMN+tHrHnV3bNQSkpKXh6ejJixAh69+7Nvffea9Wy+vbtCxhNMXFxcfj6+uLr64unpycxMTH4+Pgwbdo0wsPDcXFx4fz581y8ePGm5YSGhhIaGpqxA46Li+PYsWPExsYyYMAAvL29b1hfTh566CFWrlxJSEgI69evz0gCERERdOnShYoVKwLw8MMPEx4eDnDD9EGDBvHHH38ARsnszNdDrl69SmxsrFXbxVpmJgENpNc9KA1ccMRKTBn+7/hxUkY/x451sWx1GUSKmzvtx7Xj7hfvxrO0lBsWIjtubm7s2rWL9evXs3LlShYvXpxR9jk36U0nLi4uGf9Pf5ySksKKFSuIjIxkz549lChRgtq1a99Qljqd1pqpU6fy9NNP3zB90aJFNt2k2adPHyZNmoS/vz+lSv1b2uVWSlenpaWxY8eObMdJsBczD0fHAfOUUmeB+cDU7GZSSj2llNqtlNodGRlp80rSh//LzGHD/12/jn55Boca3c+SsIZsIJA6997BM4fHEDQvSBKAELmIi4vjypUr9OrVi0WLFmWcKfj6+ubr6PfKlSv4+flRokQJNm7cyOnTp7Ndbo8ePfjoo4+Ii4sD4Pz581y6dIm7776b1atXc/36dWJjY/nf//6X6/q8vLyYO3cu06dPv2F627Zt2bx5M1FRUaSmpvLll1/SuXNn2rZty6ZNm4iOjiY5OZlvvvkm4zVBQUEsXrw447G1F9Nt4dAzAaVUGFA5m6emA92A8Vrrb5VSDwL/BQKzzqi1XgYsA6OKqK0xpA//N2nVfhJT0hw3/N/PP/P3Uy8RfK4Jp7kPv9vLMnTxvdTtVte+6xGiGMh6TaBnz56MHTuWfv36kZCQgNaahQsXAkbzypNPPsnbb7/NqlWrbF7Xww8/TJ8+ffD396dFixbcfvvtgDHg/V133UWTJk245557mDdvHocPH6Z9+/YAlCxZks8//5xWrVoxaNAgWrRoQa1atejUqVOe63zooYdumlalShVef/11AgIC0FrTq1cv+vXrBxjdS9u3b0+VKlVo1apVxpCWb7/9Ns888wzNmjUjJSWFu+++m/fee8/mbZAb00pJK6WuAGW01loZ50JXtNa5lsW81VLS6b2DHDL839mzXBv5PBt+SmAvrfAqVYKuc3vQakQrXNyk3V8UPlJKungrSqWkLwCdgU1AV+CYo1bkkOH/kpJInb+QnTN+Jjy5A8kuHrR9xp/OM7viVdZx7XdCCGFPZiaBJ4G3lFJuQALwlImx2ERv2sQfj84m9Gwj/iGA+l2rE7S0HxUaVjA7NCGEsIlpSUBrvRXIvrNtYXXxIpeenEbI/5L5k45UqO7Bwx88wG09bzM7MiGEuCVOccdwvqWmEr/gXTb+Zz17kprh4anoOasr/mM74FrC1ezohBDilkkSyEPq9l+IGLSAzefqkaia4/9wA7q81R/v8t5mhyaEEPkmSSAnly9z7JEZhP6UShRNqNvUmx5fPIpfk0pmRyaEEHYjSSArrYma8yGhM7ZzLKk25cqkMvj9ftQf2FyGdhTCTlxdXWnatCnJycm4ubkxbNgwxo0bh4tLzt2qT506xfbt2xkyZEgBRlr8SRLI5Pq2vWx+cAkRF6pRwrUa3cffQds5A3B1l3Z/4bwcUYU3c+2gS5cuMWTIEK5cuZJRtTM7p06d4osvvpAkYGeSBIC0y1fYM3AOG9enkUB1Wnb2pevKp/Gp7Gt2aEKYqiCq8Pr5+bFs2TJat27NjBkzOH36NEOHDuXatWsALF68mA4dOjBlyhQOHz5MixYtGDZsGAMGDMh2PmEb504CWnPipeWEzPmVyJRy1K6SSI8vh1C5c0OzIxOiUCioKrx169YlLS2NS5cu4efnx7p16/D09OTYsWMMHjyY3bt3M2fOHObPn88PP/wAQHx8fLbzCds4bRKIDttL6JDl/BFZjrLuHjz4egtuf6GvtPsLkUlBVuFNL2GTnJzMmDFj2LdvH66urhlllbOydj6RO6dLAgl/xxB+3yJ27kjDjZJ06+NFuy9fwM1HKnwKkVV6Fd7MicARVXj//PNPXF1d8fPzY+bMmVSqVInffvuNtLQ0PD2z/20uXLjQqvlE7pwmCaSlpvHr+M/YsOQw8WmetKgbS7fvxlCyeT2zQxOi0CqIKrzpQyiOGTMGpRRXrlyhevXquLi4sHz58oyKmllLP+c0n7CNUyQBnZbG8spTOBPlQ02vOHq82Z6qI/uZHZYQhZ63hxuBjSoRPjnArlV400tJp3cRHTp0KBMmTABg9OjR3H///XzzzTcEBATg4+MDQLNmzXBzc6N58+YMHz48x/mEbUwrJX0rbrWUNMDuXi/hVcaDxh9PRGUafUgIZyOlpIu3olRKukD5//SK2SEIIUShI6OeCCGEE5MkIIQTKkrNwMJ6t/K5ShIQwsl4enoSHR0tiaCY0VoTHR1tc1dZp7kmIIQwVK9enXPnzhEZGWl2KMLOPD09qV69uk2vkSQghJMpUaIEderUMTsMUUhIc5AQQjgxSQJCCOHEJAkIIYQTK1J3DCulIoHTZsdhgwpAlNlBmEy2gUG2g2yDdGZshyitdc/snihSSaCoUUrt1lr7mx2HmWQbGGQ7yDZIV9i2gzQHCSGEE5MkIIQQTkySgGMtMzuAQkC2gUG2g2yDdIVqO8g1ASGEcGJyJiCEEE5MkoAQQjgxSQIOppSap5Q6opTar5RarZQqY3ZMBU0pNVApdVAplaaUKjRd4wqCUqqnUuqoUuq4UmqK2fGYQSn1kVLqklLqgNmxmEUpVUMptVEpddjyWxhrdkzpJAk43jqgida6GfAHMNXkeMxwALgPCDc7kIKklHIFlgD3AI2BwUqpxuZGZYpPgGxvVHIiKcDzWutGQDvgmcLyXZAk4GBa61CtdYrl4S+AbXVeiwGt9WGt9VGz4zBBG+C41vpPrXUSsBLoZ3JMBU5rHQ78Y3YcZtJa/6W13mv5fyxwGKhmblQGSQIF63HgZ7ODEAWmGnA20+NzFJIfvjCPUqo20BLYaW4kBhlPwA6UUmFA5Wyemq61XmuZZzrGKeGKgoytoFizDZyQymaa9Ml2YkqpksC3wDit9VWz4wFJAnahtQ7M7Xml1DDgXqCbLqY3ZuS1DZzUOaBGpsfVgQsmxSJMppQqgZEAVmitvzM7nnTSHORgSqmewAtAX611vNnxiAIVAdRXStVRSrkDDwHfmxyTMIFSSgH/BQ5rrd80O57MJAk43mLAF1inlNqnlHrP7IAKmlJqgFLqHNAe+FEpFWJ2TAXB0iFgDBCCcSHwa631QXOjKnhKqS+BHUBDpdQ5pdQTZsdkgruAoUBXy35gn1Kql9lBgZSNEEIIpyZnAkII4cQkCQghhBOTJCCEEE5MkoAQQjgxSQJCCOHEJAkIIYQTkyQghBBOTJKAEEI4MUkCQtwipVRFpdRfSqmXMk1rppRKUEo9YGZsQlhL7hgWIh+UUj2A/wGdgX3AbmCX1voxUwMTwkqSBITIJ6XUIqAvsBnoBLTQWseZG5UQ1pEkIEQ+KaU8gN+A+kAHrXWhGCxECGvINQEh8q82xrgBGqhrbihC2EbOBITIB8tAITuAYxjDBc4Ammmtz5gZlxDWkiQgRD4opeYAQ4BmwBWMMaS9gACtdZqZsQlhDWkOEuIWKaU6A88Dj2qtYyxDhw4HGmGMJidEoSdnAkII4cTkTEAIIZyYJAEhhHBikgSEEMKJSRIQQggnJklACCGcmCQBIYRwYpIEhBDCiUkSEEIIJ/Z/2GRHslAWsxMAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# generate data\n",
    "np.random.seed(123)\n",
    "N = 20\n",
    "beta0 = -4\n",
    "beta1 = 2\n",
    "x = np.random.randn(N)\n",
    "e = np.random.randn(N)\n",
    "y = beta0 + beta1*x + e\n",
    "true_x = np.linspace(min(x), max(x), 100)\n",
    "true_y = beta0 + beta1*true_x\n",
    "\n",
    "# estimate model \n",
    "beta1_hat = sum((x - np.mean(x))*(y - np.mean(y)))/sum((x - np.mean(x))**2)\n",
    "beta0_hat = np.mean(y) - beta1_hat*np.mean(x)\n",
    "fit_y = beta0_hat + beta1_hat*true_x\n",
    "\n",
    "# plot\n",
    "fig, ax = plt.subplots()\n",
    "sns.scatterplot(x, y, s = 40, label = 'Data')\n",
    "sns.lineplot(true_x, true_y, color = 'red', label = 'True Model')\n",
    "sns.lineplot(true_x, fit_y, color = 'purple', label = 'Estimated Model')\n",
    "ax.set_xlabel('x', fontsize = 14)\n",
    "ax.set_title(fr\"Linear Regression for $y = {beta0} + ${beta1}$x + \\epsilon$\", fontsize = 16)\n",
    "ax.set_ylabel('y', fontsize=14, rotation=0, labelpad=10)\n",
    "ax.legend(loc = 4)\n",
    "sns.despine()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Extensions of Ordinary Linear Regression\n",
    "\n",
    "There are many important extensions to linear regression which make the model more flexible. Those include {doc}`/content/c2/s1/regularized`—which balances the bias-variance tradeoff for high-dimensional regression models—{doc}`/content/c2/s1/bayesian`—which allows for prior distributions on the coefficients—and {doc}`/content/c2/s1/GLMs`—which introduce non-linearity to regression models. These extensions are discussed in the next chapter."
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Edit Metadata",
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
